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In  this  paper,  direct  numerical  simulation  (DNS)  data  of  an  opposed-jet  hydrogen/air  diffusion  flame 
are,  in  a  postprocessing  step,  analyzed  using  the  proper  orthogonal  decomposition  (POD)  technique.  The 
aim  of  this  work  is  twofold.  The  first  goal  is  to  compute  a  small  number  of  space-dependent  empirical 
eigenfunctions,  so  that  a  low-dimensional  representation  of  the  data  generated  by  the  large  model  of  the 
discretized  partial  differential  equations  can  be  obtained  using  a  weighted  sum  of  these  few  eigenfunctions 
(POD  modes).  It  is  found  that  only  six  modes  are  needed  for  an  accurate  representation  of  the  data  in  an 
extended  range  of  inflow  velocities.  This  large  data  reduction  takes  into  account  not  only  chemical  kinetics 
bnt  also  transport  phenomena  in  a  full  two-dimensional  context  and  constitutes  the  first  step  toward  the 
construction  of  low-dimensional  dynamic  models  for  the  opposed-jet  system.  It  is  also  found  that  the  PODs 
have  very  good  interpolatory  properties.  The  second  goal  is  to  use  part  of  the  available  data  (i.e.,  partial 
measurements),  together  with  the  computed  modes,  to  estimate,  or,  in  the  terminology  of  process  control, 
to  observe,  the  “unmeasured”  quantities.  It  is  found  that  only  a  small  number  of  measurements  are  needed 
to  obtain  accurate  estimates  of  the  rest  of  the  data. 

Introduction  literature  [3])  has  been  used  in  the  construction  of 

low-dimensional  models  from  extensive  numerical 
An  accurate,  well-resolved  discretization  of  the  or  experimental  data  from  non-reactive  turbulent 
conservation  equations  for  laminar  chemically  reac-  flows  (see  Ref.  [4]  for  a  detailed  presentation  of  the 
tive  flows  reduces  the  partial  differential  equations  method  and  its  application). 

(PDEs,  which  are  infinite-dimensional  dynamical  In  this  paper,  we  employ  the  POD  technique  to 
systems)  to  systems  of  ordinary  differential  equa-  obtain  low-dimensional  representations  of  a  H2/air 
tions  (ODEs).  However,  the  size  of  the  latter  model  opposed-jet  diffusion  flame  over  a  range  of  jet-exit 
is  usually  veiy  large:  in  two  dimensions  and  with  de-  velocities.  The  statistical  analysis  of  the  simulation 
tailed  chemical  kinetics,  such  models  can  easily  results  provides  an  optimal  basis  of  orthogonal  ei- 
reach  (and  even  exceed)  O  (105)  degrees  of  freedom.  genfunctions  (modes),  in  the  sense  that  they  mini- 
While  the  descriptive  and  predictive  capabilities  of  mize  the  mean  square  error  between  the  original 
such  models  for  a  large  range  of  operating  parame-  data  and  its  low-dimensional  representation.  It  turns 
ters  can  be  exceptional,  the  computational  cost  can  out  that  the  reactive-flow  system  with  more  than 
also  be  prohibitive.  26,000  degrees  of  freedom  can  be  approximated  by 

The  method  of  proper  orthogonal  decomposition  the  weighted  sum  of  only  six  inodes.  This  in  turn 
(POD)  (also  known  as  the  method  of  empirical  or-  means  that  one  only  needs  the  six  coefficients 
thogonal  eigenfunctions  [EOFs]  in  the  meteorology  (weights)  of  the  modes  to  get  a  good  description  of 
and  oceanography  literature  [I],  as  Karhunen-Loeve  the  system.  We  find  that  this  low-dimensional  rep- 
[KL]  expansion  in  pattern  recognition  [2],  and  as  resentation  approximates  well  not  only  the  direct  nu- 
principal  component  analysis  [PCA]  in  the  statistical  merical  simulation  (DNS)  data  used  to  obtain  the 
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modes  but  also  data  for  intermediate  jet-exit  veloc¬ 
ities  within  the  range  of  velocities  covered  by  the 
data. 

In  principle,  one  can  substitute  the  original  vari¬ 
ables  by  their  weighted  sum  expansion  in  the  PDEs 
and  perform  a  Galerkin  procedure  to  obtain  a  system 
of,  here,  only  six  ODEs  for  the  time  derivatives  of 
the  coefficients  in  the  expansion  (a  POD-Galerkin 
procedure).  Time  integration  of  this  system  could 
then  provide  a  reasonably  accurate  approximation  to 
the  full  (more  than  26,000  ODEs)  system.  For  sim¬ 
ple  forms  of  the  nonlinear  terms  in  the  original 
PDEs  (e.g.,  quadratic),  the  low-dimensional  model 
can  be  constructed  in  explicit  form  [5,6].  However, 
for  higher-order  and/ or  non-polynomial  nonlineari¬ 
ties  (as  is  the  case  for  reactive  systems  due  to  the 
exponential  temperature  dependence  of  the  reaction 
rate  constants),  the  evaluation  of  the  right-hand  side 
and  the  Jacobian  for  POD-Galerkin  models  is  com¬ 
putationally  expensive. 

Plere  we  are  interested  in  another  problem:  the 
construction  of  POD-based  observers,  so  that  the 
computed  POD  modes  and  a  set  of  partial  measure¬ 
ments  can  be  used  to  estimate  the  complete  field  of 
unmeasured  quantities.  Typically  in  experiments, 
measurements  are  (or  can  only  be)  made  for  a  few 
variables  either  at  selected  positions  (e.g.,  velocity 
with  laser  Doppler  anemometry,  temperature  with 
coherent  anti-Raman  Stokes  or  a  thermocouple)  or, 
in  the  best  case,  over  two-dimensional  domains  (e.g., 
velocity  with  particle  image  velocimetry  [PIV],  OH 
concentration  with  planar  laser-induced  fluores¬ 
cence  [PLIF],  temperature  with  Rayleigh  scatter¬ 
ing).  In  particular,  opposed-jet  diffusion  flames  are 
commonly  investigated  experimentally  by  measuring 
the  velocity  and  temperature  at  points  along  the  axis 
of  symmetry.  In  more  recent  experiments  investi¬ 
gating  transient  processes  such  as  flame-vortex  in¬ 
teractions  in  the  same  setup,  two-dimensional  mea¬ 
surement  techniques  such  as  PLIF  were  employed. 
Recently,  Podvin  and  Lumley  [7]  used  POD  modes 
for  the  velocity  in  the  wall  region  of  a  turbulent 
boundary  layer  to  reconstruct  the  flow  in  the  wall 
region  only  from  wall-shear  measurements.  We  in¬ 
vestigate  this  idea  here  to  find  how  well  partial  mea¬ 
surements  of  a  single  variable  (either  in  the  full  do¬ 
main  or  along  the  axis  of  symmetry)  can  be  used 
together  with  the  available  POD  modes  to  estimate 
the  unmeasured  variables;  the  results  are  remarkably 
good. 

Methodology  and  Results 

Simulation  of  the  Opposed-Jet  Hz/Air  Burner 

The  DNS  data  were  obtained  with  a  spectral-ele¬ 
ment  code  that  solves  the  low  Mach  number  form 
of  the  conservation  equations  for  momentum,  spe¬ 


cies,  and  energy  of  compressible  chemically  reactive 
flows  [8],  The  problem  solved  was  that  of  a  20%  H2 
(by  mole)  in  N2  and  air  diffusion  flame,  stabilized  on 
an  opposed-jet  burner  in  a  two-dimensional,  axisym- 
metric  geometiy.  The  geometry,  the  boundary  con¬ 
ditions,  and  the  grid  employed  in  the  simulations  are 
shown  in  Fig.  1.  For  fixed  geometry  inflow-stream 
compositions  and  temperature,  and  for  equal  veloc¬ 
ity  magnitudes  of  the  plug-flow  profiles  at  the  jet 
exits,  the  behavior  of  the  system  is  parameterized 
only  by  the  Reynolds  number  defined  at  the  air  side, 
Re  =  wa jr  *  D/v^.  Extending  the  work  presented  in 
Ref.  [9],  simulations  were  performed  in  the  range 
from  low-  to  high-flow-rate  extinction  of  the  flame, 
corresponding  to  Reynolds  number  values  from  Re 
=  10  to  Re  =  1,800.  For  each  Reynolds  number, 
the  simulation  starts  from  an  appropriate  initial  con¬ 
dition,  and  the  whole  transient  response,  until  the 
steady-state  flame  is  established,  is  computed  by 
time-accurate  integration.  Resolution  is  adjusted  to 
accommodate  the  flame  thickness  by  varying  the  or¬ 
der  of  the  interpolating  polynomial  in  each  element. 
For  the  low-  to  medium-strain-rate  simulations  (Re 
—  20,  50,  100,  200,  and  400),  it  was  found  that 
fourth-order  Legendre  polynomials  (5  X  5  elemen¬ 
ted  resolution)  in  both  the  radial  and  axial  direction 
of  the  120  spectral-element  mesh  provide  an  accu¬ 
rate  discretization.  The  detailed  chemical  kinetics  of 
Ref.  [10]  (9  species,  19  reversible  elementary  reac¬ 
tions)  and  multicomponent  transport  based  on  the 
mixture-fraction  formulation  were  used  in  the  sim¬ 
ulations. 

Data  Reduction 

Forty-nine  snapshots  (two-dimensional  fields  at 
steady  state  and  along  the  transients)  of  velocity, 
pressure,  temperature,  and  species  mass  fractions  (a 
total  of  13  variables  each)  were  extracted  from  the 
simulations  and  interpolated  on  a  uniform  101  X 


axis 

I 


Fig.  1.  Spectral  element  skeleton,  coordinate  system, 
and  boundary  conditions  used  for  the  opposed-jet  DNS 
(uniform  inflow  velocity  profiles  for  both  fuel  and  air). 
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201  grid.  To  account  for  the  different  units  and  mag¬ 
nitudes  of  the  variables,  they  were  all  normalized 
between  0  and  1  and  then  arranged  in  a  matrix, 
V  =  [u],  «2,  .  .  .  ,  w^],  of  M  =  49  one-dimen¬ 
sional  column  vectors  of  length  101  X  201  X  13. 
The  latter  was  subsequently  processed  with  the  code 
KLTOOL  [11],  which  found  the  orthogonal  eigen¬ 
vectors  (the  POD  modes)  using  the  method  of  snap¬ 
shots  [12]  as  follows: 

1.  Compute  the  mean  vector,  u,  with 


2.  Compute  the  vectors  with  zero  mean  (variations), 
ii;  =  u'j  —  u,  and  arrange  them  in  a  matrix,  U, 
as  before. 

3.  Compute  the  M  X  M  covariance  matrix, 

c  =  i  (uTu). 

M 


4.  Compute  the  non-negative  eigenvalues,  1;,  i  = 
1,  .  .  .  ,  M,  in  the  order  of  decreasing  magnitude 
and  the  corresponding  eigenvectors,  </>;,  of  the  co- 
variance  matrix,  C. 

The  orthogonal  eigenfunctions  (modes)  can  then  be 
written  as 

M 

y^k  &ik  D,' 

i  =  1 

where  is  the  ith  component  of  the  fcth  eigenvec¬ 
tor.  As  an  example,  the  mean  field  (dimensional)  and 
the  six  modes  (dimensionless)  of  temperature  are 
shown  in  Fig.  2.  The  computed  modes  describe  the 
variations  of  the  data  from  the  mean.  As  Re  is  in¬ 
creased,  the  flame  becomes  thinner  and  longer,  as 
reflected  by  the  isocontours  of  the  temperature 
modes  that  show  large  variations  on  both  sides  and 
close  to  the  tip  of  the  mean  field.  By  projecting  the 
database  vectors  onto  the  modes  y/^x),  k  —  1,  .  .  ., 
N  <  M,  the  coefficient,  tfj.,  that  provide  the  optimal 
approximation  to  the  full  data, 


T  (mean) 


0  0.5  1  1.5  2 


r  [cm] 

T  (mode  1  )  T  (mode  2  ) 


Fig.  2.  Mean  field  and  first  six  POD  modes  for  temperature.  The  mean  field  is  dimensional,  while  the  modes  are 
dimensionless. 


78 


TURBULENT  COMBUSTION — Direct  Numerical  Simulations  and  Large  Eddy  Simulations 


Se+04 
7e+04 
6e+04 
g  5c+04 
S  4e+04 
3e+04 
2e+04 
lc+04 
0e+00 


Fig.  3.  Histogram  of  the  nine  largest  eigenvalues  (left 
axis)  and  percentage  of  the  energy  captured  by  the  corre¬ 
sponding  modes  (line  with  filled  circles,  right  axis). 


N 

u'  =  U  +  2  k'k^VkM 

k  =  l 

are  obtained  as 


ak  ~ 


» ■  Mb. 

Vk '  'Vk 


The  sum  of  the  eigenvalues  of  the  covariance  matrix 
is  defined  as  the  “energy”  of  the  data.  A  commonly 
used  criterion  for  selecting  the  number  of  modes  to 
keep  in  a  low-dimensional  representation  is  that  of 
capturing  a  certain  percentage  of  the  energy  of  the 
data.  As  can  be  seen  in  Fig.  3,  only  N  —  6  and  N  = 
9  modes  are  needed  to  capture  2)1  =  1  An/S„=1  Xn  — 
99%  and  99.5%  of  the  energy  of  the  DNS  data,  re¬ 
spectively.  The  rest  of  the  modes  are  considered 
comparatively  meaningless  and  are  ignored  from  the 
low-dimensional  representation.  Indeed,  velocities, 
temperature,  and  species  mass  fractions  (with  the 
exception  of  H02  and  H202,  which  can  be  as  much 
as  25%  off  in  some  regions)  of  all  the  original  data 
fields  are  reconstructed  with  less  than  13%  relative 
error  over  the  whole  two-dimensional  domain,  with 
both  six  and  nine  modes.  The  additional  three  modes 
offer  marginal  improvement  in  the  accuracy  of  the 
approximation  of  velocity  and  mass  fraction  of  major 
species  but  have  a  stronger  effect  on  the  mass  frac¬ 
tions  of  minor  species. 

Accurate  approximations  were  obtained  not  only 
for  the  DNS  data  used  to  compute  the  modes  but 
also  for  data  at  intermediate  values  of  the  Reynolds 
number  (i.e.,  for  Re  =  20  -f  400).  Fig.  4  compares 
the  DNS  fields  of  temperature  and  OH  mass  frac¬ 
tion  with  those  obtained  with  six  modes  for  the 
steady-state  solution  at  Re  =  300  (which  was  not 
used  to  obtain  the  POD  modes).  In  this  case,  in¬ 
creasing  the  number  of  modes  to  nine  did  not  re¬ 
duce  the  relative  error. 
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Fig.  4.  Re  =  300:  comparison  of  the  DNS  two-dimen¬ 
sional  fields  of  temperature  (top)  and  OH  mass  fractions 
(bottom)  and  their  approximation  via  six  POD  modes.  The 
distribution  of  the  relative  error  is  also  plotted. 


We  can  conclude  that  a  linear  combination  of  only 
six  modes  provides  an  accurate  low-dimensional  rep¬ 
resentation  of  the  data  generated  by  the  full  (26,000 
degrees  of  freedom)  model  for  the  range  of  Re  val¬ 
ues  used  to  compute  the  POD  modes.  In  principle, 
one  can  substitute  this  approximation  back  into  the 
conservation  equation  and  derive  the  equations  of 
change  for  the  coefficients,  a ^  (see,  for  example, 
Refs.  [5,6]  for  non-reactive  flows),  which  will  effec¬ 
tively  reduce  the  original  26,000  ODEs  to  only  6 
ODEs.  In  the  next  subsection,  however,  we  consider 
the  problem  of  using  the  available  POD  modes  to¬ 
gether  with  measurements  of  some  of  the  variables 
(e.g.,  temperature)  to  obtain  an  estimate  of  the  un¬ 
measured  variables. 
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Fig.  5.  Re  =  300:  comparison  of  the  DNS  OH  field  with 
its  estimate  using  six  (top)  and  nine  (bottom)  POD  modes 
together  with  temperature  two-dimensional  measure¬ 
ments.  The  distribution  of  the  relative  error  is  also  plotted. 

POD-Based  Observers 

Since  only  a  small  number  of  POD  mode  can  pro¬ 
vide  a  good  approximation  to  the  DNS  model,  in 
principle,  only  the  same  small  number  of  measure¬ 
ments  at  the  proper  locations  (i.e.,  with  non-trivial 
components  in  the  span  of  the  POD  modes)  is  suf¬ 
ficient  to  evaluate  the  coefficients,  a^.  from  which 
the  entire  field  can  be  approximated.  In  fact,  one 
can  even  find  the  optimal  measurement  points,  an 
issue  we  are  currently  working  on  but  which  will  not 
be  considered  here.  Therefore,  a  few  measurements 
of  easily  measured  quantities  (e.g.,  temperature  or 
velocity)  at  selected  locations,  together  with  the 
POD  modes,  can  be  used  to  effectively  obtain  in¬ 
formation  on  other  quantities  in  the  entire  domain 


which  cannot  be  easily  measured,  such  as  radical 
species  concentrations.  Typically,  measurements  are 
performed  at  many  more  points  than  the  small  num¬ 
ber  of  POD  coefficients,  resulting  in  an  overdeter¬ 
mined  system  of  linear  equations.  In  this  case,  the 
estimate  of  the  POD  coefficients  are  found  as  a  least- 
square  solution  to  the  problem  using  singular  value 
decomposition  [13], 

This  procedure  was  applied  to  our  data,  substitut¬ 
ing  the  experiment  with  direct  simulation.  Two  cases 
commonly  used  in  the  experimental  study  of  op¬ 
posed-jet  burners  were  considered:  measurement 
(1)  of  two-dimensional  fields  and  (2)  along  the  axis 
of  symmetry.  Fig.  5  compares  the  computed  OH 
mass  fraction  fields  with  those  obtained  using  six 
modes  (top  half)  and  nine  modes  (bottom  half)  and 
the  temperature  two-dimensional  field  measure¬ 
ment  (i.e.,  at  101  X  201  points)  to  estimate  the  co¬ 
efficients,  Ojj,  for  the  Re  =  300  case.  It  should  be 
stressed  that  the  Re  =  300  data  were  not  used  to 
compute  the  POD  modes.  Again,  with  the  exception 
of  HOa  and  II202,  the  “unmeasured”  quantities  can 
be  estimated  with  less  than  13%  maximum  relative 
error.  In  this  case,  increasing  the  number  of  modes 
to  nine  offered  only  marginal  improvement.  Good 
results  were  obtained  even  when  temperature  was 
measured  at  only  101  points  along  the  axis  of  sym¬ 
metry,  as  shown  in  Fig.  6  for  the  Re  =  300  case. 
The  behavior  and  values  of  the  relative  errors  are 
practically  unaffected  in  a  large  region  around  the 
axis.  This  indicates  that  it  is  indeed  possible  to  use 
only  a  handful  of  carefully  selected  measurement 
points  (together  with  the  POD  modes)  to  estimate 
the  remaining  quantities  over  the  two-dimensional 
domain.  Higher  relative  errors  are  observed  in  the 
low  OH  mass  fraction  regions  away  from  the  axis, 
where,  in  the  case  of  axial  data  only,  the  observer  is 
not  expected  to  perform  well  in  the  first  place.  When 
the  number  of  modes  is  increased  from  six  (Fig.  6, 
left  half)  to  nine  (Fig.  6,  right  half)  close  to  the  axis, 
the  relative  error  decreases.  The  relative  error  is  al¬ 
most  doubled,  however,  in  a  region  one  diameter 
away  from  the  axis,  where,  again,  the  OH  mass  frac¬ 
tion  is  very  low.  In  passing,  we  would  like  to  com¬ 
ment  on  an  issue  raised  by  a  reviewer  with  respect 
to  the  optimal  number  of  modes  for  the  observers. 
Different  criteria  have  been  proposed  [1],  with  the 
energy  criterion  used  here  being  the  most  popular 
one.  Once  the  energy  threshold  is  chosen,  the  opti¬ 
mal  representation  is  fixed.  The  quality  of  the  ob¬ 
server  estimates  is  then  solely  determined  by  the 
principal  angles  [13]  between  the  POD  and  the  mea¬ 
surement  subspace.  It  should  also  be  noted  that  the 
maximum  relative  error  (L„  norm)  is  just  one  way  to 
quantify  the  accuracy  of  the  estimations;  any  error 
norm  can  be  selected  (in  fact,  the  L2  norm  would 
give  lower  error  values).  We  would  like  to  add  at  this 
point,  as  Prof.  Gouldin  pointed  out  during  the  sym¬ 
posium  to  one  of  the  authors  (C.E.F.),  that  the  Kar- 
hunen-Loeve  basis  was  exploited  in  a  similar  way  in 
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Fig.  6.  Re  =  300:  comparison  of  the  DNS  OH  field  with 
its  estimate  using  six  (top)  and  nine  (bottom)  POD  modes 
together  with  temperature  measurements  along  the  axis  of 
symmetry.  The  distribution  of  the  relative  error  is  also  plot¬ 
ted. 


the  evaluation  of  a  tomographic  inversion  method  as 
a  combustion  diagnostic  tool  [14], 

Conclusions  and  F urther  Work 

The  method  of  POD  was  applied  to  DNS  data  for 
an  opposed-jet  H2/air  burner.  It  was  found  that  only 
six  modes  capture  99%  of  the  total  “energy”  of  the 
original  system  of  26,000  degrees  of  freedom.  In 
contrast  to  the  traditional  reduction  techniques  ap¬ 
plied  to  chemically  reactive  systems,  this  method 
takes  into  account  not  only  chemical  kinetics  but  also 
transport  and  transient  phenomena  in  a  full  two- 
dimensional  context.  The  results  suggest  that  the 


problem  investigated  here  (26,000  ODEs)  can  be 
described  effectively  with  only  6  ODEs.  However, 
the  construction  of  a  low-dimensional  dynamic 
model  based  on  the  POD  modes  is  complicated  by 
the  highly  nonlinear  form  of  the  reaction  rate  terms 
in  the  conservation  equations.  It  should  be  possible 
to  overcome  this  problem  either  by  moving  between 
the  (small)  POD  space  and  the  original  (large)  space 
or  by  fitting  a  functional  form  to  the  right-hand  side 
of  the  ODEs  for  the  POD  coefficients.  We  plan  to 
explore  both  approaches,  which  will  allow  for  the 
true  model  reduction  to  only  six  ODEs,  in  the  near 
future. 

In  addition  to  data  reduction,  the  POD  modes 
were  used  in  a  way  which  can  be  of  interest  in  ex¬ 
periments.  Since  only  a  few  modes  are  needed  to 
provide  a  good  representation  of  the  full  system, 
only  a  handful  of  measurements  at  selected  points 
can,  in  principle,  be  used  to  obtain  estimates  of  the 
coefficients,  ah  and  therefore  of  the  unmeasured 
quantities.  The  problem  of  selecting  the  measure¬ 
ment  points  was  not  addressed  in  this  paper.  Instead, 
the  good  performance  of  the  proposed  methodology 
was  shown  for  the  case  where  the  number  of  mea¬ 
surements  exceeds  the  number  of  modes.  The  co¬ 
efficients  in  this  case  are  obtained  as  a  least-square 
solution  to  the  corresponding  overdetermined  sys¬ 
tem  of  linear  equations.  We  are  currently  working 
on  an  algorithm  for  the  selection  of  the  minimal 
number  of  sensors  and  their  location  for  the  accurate 
estimation  of  the  complete  field.  Since  well-resolved 
DNS  data  were  used  to  obtain  the  POD  modes,  we 
did  not  address  here  the  effect  of  noise  in  the  initial 
data.  This  will  be  particularly  important  for  experi¬ 
mental  data  processing,  and  we  plan  to  examine  this 
in  the  future. 
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